This script is similar to the ek_irrad_model but uses relative Hsat for analysis The amount of time in Hsat was divided by the total day length, which is defined by the amount of time when umols photons m-2 s-1 were greater than or equal to 1.
#load the various libraries
library(lme4)
library(lmerTest)
library(effects)
library(car)
library(MuMIn)
library (dplyr)
library(emmeans)
library(DHARMa)
library(performance)
library(patchwork)
library(rstatix)
#for plots and tables
library(ggplot2)
library(ggpubr)
library(forcats)
library(RColorBrewer)
library(tidyverse)
library(sjPlot)
library(sjmisc)
library(mmtable2)
library(gt)
library(purrr)
library(stringr)
library(tidyr)
Load the dataset
ek_irrad_data <- read.csv("input_data/mean_supersaturation_by_period.csv")
Assign various variables as factors
ek_irrad_data$Run <- as.factor(ek_irrad_data$Run)
ek_irrad_data$Temperature <- as.factor(ek_irrad_data$Temp...C.)
ek_irrad_data$Treatment <- as.factor(as.character(ek_irrad_data$Treatment))
#ek_irrad_data$deltaNPQ <- as.factor(ek_irrad_data$deltaNPQ)
ULVA____________________________________________________________
Subset the data by species and make new column of the treatments for the plots
ulva <- subset(ek_irrad_data, Species == "ul" & Treatment != 2.5)
ulva$treatment_graph[ulva$Treatment == 0] <- "1) 35ppt/0.5umol"
ulva$treatment_graph[ulva$Treatment == 1] <- "2) 35ppt/14umol"
ulva$treatment_graph[ulva$Treatment == 2] <- "3) 28ppt/27umol"
ulva$treatment_graph[ulva$Treatment == 3] <- "5) 18ppt/53umol"
ulva$treatment_graph[ulva$Treatment == 4] <- "6) 11ppt/80umol"
ulva$treatment_graph[ulva$Treatment == 2.5] <- "4) 28ppt/53umol"
Make a histogram for ulva
hist(ulva$supersat_rel, main = paste("Ulva lactuca"), col = "olivedrab3", labels = TRUE)
ulva %>% ggplot(aes(supersat_avg)) +
geom_histogram(binwidth=5, fill = "#5BB300", color = "black", size = 0.25, alpha = 0.85) +
theme_bw()
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
Run model without interaction between the treatments and temperature Take RLC.Order out of the random effects because causing problems of singularity. R2 is same with or without (+ (1 | RLC.Order))
rel_hsat_model_ulva <- lmer(formula = supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = ulva)
Make residual plots of the data
hist(resid(rel_hsat_model_ulva))
plot(resid(rel_hsat_model_ulva) ~ fitted(rel_hsat_model_ulva))
qqnorm(resid(rel_hsat_model_ulva))
qqline(resid(rel_hsat_model_ulva))
Check the performance of the model, get R2, plot the effects of model and make table
performance::check_model(rel_hsat_model_ulva)
## Variable `Component` is not in your data frame :/
r.squaredGLMM(rel_hsat_model_ulva)
## Warning: 'r.squaredGLMM' now calculates a revised statistic. See the help page.
## R2m R2c
## [1,] 0.1917698 0.822725
summary(rel_hsat_model_ulva)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) +
## (1 | RLC.Order)
## Data: ulva
##
## REML criterion at convergence: 1443.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.8397 -0.4267 0.0763 0.5141 2.6191
##
## Random effects:
## Groups Name Variance Std.Dev.
## Plant.ID (Intercept) 14.3670 3.790
## Run (Intercept) 39.1517 6.257
## RLC.Order (Intercept) 0.9662 0.983
## Residual 15.3082 3.913
## Number of obs: 240, groups: Plant.ID, 96; Run, 7; RLC.Order, 6
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 74.4406 4.5814 5.2312 16.248 1.12e-05 ***
## Treatment1 -5.8265 5.3559 4.9914 -1.088 0.3264
## Treatment2 -5.2599 5.3559 4.9914 -0.982 0.3712
## Treatment3 -10.2494 5.3559 4.9914 -1.914 0.1139
## Treatment4 -11.1703 5.3559 4.9914 -2.086 0.0915 .
## Temperature27 0.9231 1.4722 5.3816 0.627 0.5563
## Temperature30 1.6739 1.3225 24.8618 1.266 0.2173
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Trtmn1 Trtmn2 Trtmn3 Trtmn4 Tmpr27
## Treatment1 -0.823
## Treatment2 -0.823 0.989
## Treatment3 -0.823 0.989 0.989
## Treatment4 -0.823 0.989 0.989 0.989
## Temperatr27 -0.151 -0.001 -0.001 -0.001 -0.001
## Temperatr30 -0.145 0.000 0.000 0.000 0.000 0.460
plot(allEffects(rel_hsat_model_ulva))
tab_model(rel_hsat_model_ulva, show.intercept = TRUE, show.se = TRUE, show.stat = TRUE, show.df = TRUE, show.zeroinf = TRUE)
| supersat_rel | ||||||
|---|---|---|---|---|---|---|
| Predictors | Estimates | std. Error | CI | Statistic | p | df |
| (Intercept) | 74.44 | 4.58 | 65.41 – 83.47 | 16.25 | <0.001 | 229.00 |
| Treatment [1] | -5.83 | 5.36 | -16.38 – 4.73 | -1.09 | 0.278 | 229.00 |
| Treatment [2] | -5.26 | 5.36 | -15.81 – 5.29 | -0.98 | 0.327 | 229.00 |
| Treatment [3] | -10.25 | 5.36 | -20.80 – 0.30 | -1.91 | 0.057 | 229.00 |
| Treatment [4] | -11.17 | 5.36 | -21.72 – -0.62 | -2.09 | 0.038 | 229.00 |
| Temperature [27] | 0.92 | 1.47 | -1.98 – 3.82 | 0.63 | 0.531 | 229.00 |
| Temperature [30] | 1.67 | 1.32 | -0.93 – 4.28 | 1.27 | 0.207 | 229.00 |
| Random Effects | ||||||
| σ2 | 15.31 | |||||
| τ00 Plant.ID | 14.37 | |||||
| τ00 Run | 39.15 | |||||
| τ00 RLC.Order | 0.97 | |||||
| ICC | 0.78 | |||||
| N Run | 7 | |||||
| N Plant.ID | 96 | |||||
| N RLC.Order | 6 | |||||
| Observations | 240 | |||||
| Marginal R2 / Conditional R2 | 0.192 / 0.823 | |||||
Construct null model to perform likelihood ratio test REML must be FALSE
ulva_relhsat_treatment_null <- lmer(formula = supersat_rel ~ Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = ulva, REML = FALSE)
ulva_relhsat_model2 <- lmer(formula = supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = ulva, REML = FALSE)
anova(ulva_relhsat_treatment_null, ulva_relhsat_model2)
## Data: ulva
## Models:
## ulva_relhsat_treatment_null: supersat_rel ~ Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order)
## ulva_relhsat_model2: supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order)
## npar AIC BIC logLik deviance Chisq Df
## ulva_relhsat_treatment_null 7 1544.8 1569.1 -765.39 1530.8
## ulva_relhsat_model2 11 1481.6 1519.9 -729.81 1459.6 71.146 4
## Pr(>Chisq)
## ulva_relhsat_treatment_null
## ulva_relhsat_model2 1.3e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ulva_relhsat_temperature_null <- lmer(formula = supersat_rel ~ Treatment + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = ulva, REML = FALSE)
ulva_relhsat_model3 <- lmer(formula = supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = ulva, REML = FALSE)
anova(ulva_relhsat_temperature_null, ulva_relhsat_model3)
## Data: ulva
## Models:
## ulva_relhsat_temperature_null: supersat_rel ~ Treatment + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order)
## ulva_relhsat_model3: supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order)
## npar AIC BIC logLik deviance Chisq Df
## ulva_relhsat_temperature_null 9 1479.0 1510.4 -730.53 1461.0
## ulva_relhsat_model3 11 1481.6 1519.9 -729.81 1459.6 1.4286 2
## Pr(>Chisq)
## ulva_relhsat_temperature_null
## ulva_relhsat_model3 0.4895
Plots
ulva %>% ggplot(aes(treatment_graph, supersat_rel)) +
geom_boxplot(size=0.5) +
geom_point(alpha = 0.75, size = 3, aes(color = Temperature), show.legend = TRUE) +
labs(x="Treatment (salinty/nitrate)", y= "Relative Hsat (%)", title= "C", subtitle = "Ulva lactuca") +
scale_x_discrete(labels = c("35ppt/0.5umolN", "35ppt/14umolN", "28ppt/27umolN", "18ppt/53umolN", "11ppt/80umolN")) +
ylim(20, 100) + stat_mean() +
geom_hline(yintercept=60, color = "red", size = 0.5, alpha = 0.5) +
scale_color_manual(values = c("#295102", "#7CB950", "#BDE269")) +
theme_bw() +
theme(legend.position = c(0.90,0.90), plot.title = element_text(face = "bold", vjust = -15, hjust = 0.05),
plot.subtitle = element_text(face = "italic", size = 14, vjust = -18, hjust = 0.05))
Summarize the means for relative Hsat
ulva %>% group_by(Treatment) %>% summarise_at(vars(supersat_rel), list(mean = mean))
## # A tibble: 5 × 2
## Treatment mean
## <fct> <dbl>
## 1 0 75.3
## 2 1 70.1
## 3 2 70.7
## 4 3 65.7
## 5 4 64.7
ulva %>% group_by(Run) %>% summarise_at(vars(supersat_rel), list(mean = mean))
## # A tibble: 7 × 2
## Run mean
## <fct> <dbl>
## 1 1 72.3
## 2 2 73.8
## 3 3 69.2
## 4 3.5 60.0
## 5 4 60.4
## 6 6 71.2
## 7 6.5 79.4
ulva %>% group_by(Run) %>% summarise_at(vars(day_length_avg), list(mean = mean))
## # A tibble: 7 × 2
## Run mean
## <fct> <dbl>
## 1 1 654.
## 2 2 643.
## 3 3 626.
## 4 3.5 621.
## 5 4 611.
## 6 6 626.
## 7 6.5 633.
Add growth rate to this dataset
growth_rate <- read.csv("/Users/Angela/src/work/limu/algal_growth_photosynthesis/data_input/all_runs_growth_011723.csv")
growth_rate$Species <- as.factor(growth_rate$Species)
growth_rate$treatment <- as.factor(growth_rate$treatment)
Make a new column for weight change (difference final from initial)
growth_rate$growth_rate_percent <- (growth_rate$final.weight - growth_rate$Initial.weight) / growth_rate$Initial.weight * 100
Subset by species and add to current dataset
gr_ulva <- subset(growth_rate, Species == "Ul" & treatment != 2.5)
ulva$growth_rate <- round((gr_ulva$final.weight - gr_ulva$Initial.weight) / gr_ulva$Initial.weight * 100, digits = 2)
Plot a regression between the photosynthetic independent variables of interest and growth rate
ulva_growth_rel_hsat_graph <- ggplot(ulva, aes(x=supersat_rel, y=growth_rate)) +
geom_point(alpha = 0.5, size = 3, aes(color = Treatment), show.legend = TRUE) +
geom_smooth(method = "lm", col = "black") + theme_bw() +
labs(title = "A", subtitle = "Ulva lactuca", x = "Mean Rel. Hsat Time (%)",
y = "growth rate (%)") + stat_regline_equation(label.x = 75, label.y = 225) +
stat_cor(label.x = 75, label.y = 215) +
theme(plot.title = element_text(face = "bold", vjust = -15, hjust = 0.05),
plot.subtitle = element_text(face = "italic", vjust = -20, hjust = 0.05))
ulva_growth_rel_hsat_graph
## `geom_smooth()` using formula = 'y ~ x'
HYPNEA_____________________________________________________________________________________________________________ Subset for Hypnea
hypnea <- subset(ek_irrad_data, Species == "hm" & day1_rlc_time != "11:34:05")
hypnea$treatment_graph[hypnea$Treatment == 0] <- "1) 35ppt/0.5umol"
hypnea$treatment_graph[hypnea$Treatment == 1] <- "2) 35ppt/14umol"
hypnea$treatment_graph[hypnea$Treatment == 2] <- "3) 28ppt/27umol"
hypnea$treatment_graph[hypnea$Treatment == 3] <- "5) 18ppt/53umol"
hypnea$treatment_graph[hypnea$Treatment == 4] <- "6) 11ppt/80umol"
hypnea$treatment_graph[hypnea$Treatment == 2.5] <- "4) 28ppt/53umol"
Make a histogram of the data for ulva
hist(hypnea$supersat_rel, main = paste("Hypnea musciformis"), col = "maroon", labels = TRUE)
hypnea %>% ggplot(aes(supersat_rel)) +
geom_histogram(binwidth=5, fill = "#a8325e", color = "black", size = 0.25, alpha = 0.85) +
theme_bw()
Run model
rel_hsat_model_hypnea <- lmer(formula = supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = hypnea)
Plot residuals
hist(resid(rel_hsat_model_hypnea))
plot(resid(rel_hsat_model_hypnea) ~ fitted(rel_hsat_model_hypnea))
qqnorm(resid(rel_hsat_model_hypnea))
qqline(resid(rel_hsat_model_hypnea))
Check the performance of the model, get R2, make table and plot effects
performance::check_model(rel_hsat_model_hypnea)
## Variable `Component` is not in your data frame :/
r.squaredGLMM(rel_hsat_model_hypnea)
## R2m R2c
## [1,] 0.2345395 0.7173404
summary(rel_hsat_model_hypnea)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) +
## (1 | RLC.Order)
## Data: hypnea
##
## REML criterion at convergence: 1824.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3952 -0.4447 0.0417 0.5604 2.3878
##
## Random effects:
## Groups Name Variance Std.Dev.
## Plant.ID (Intercept) 15.7258 3.9656
## Run (Intercept) 24.7192 4.9718
## RLC.Order (Intercept) 0.8515 0.9228
## Residual 24.1774 4.9171
## Number of obs: 286, groups: Plant.ID, 144; Run, 8; RLC.Order, 6
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 65.012 3.732 5.286 17.419 7.14e-06 ***
## Treatment1 -7.659 4.361 5.043 -1.756 0.138852
## Treatment2 -5.662 4.361 5.043 -1.298 0.250378
## Treatment2.5 -10.804 6.224 4.561 -1.736 0.148734
## Treatment3 -7.924 4.361 5.043 -1.817 0.128380
## Treatment4 -10.498 4.363 5.056 -2.406 0.060613 .
## Temperature27 4.541 1.410 4.683 3.221 0.025750 *
## Temperature30 6.123 1.310 13.044 4.674 0.000432 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Trtmn1 Trtmn2 Trt2.5 Trtmn3 Trtmn4 Tmpr27
## Treatment1 -0.811
## Treatment2 -0.811 0.974
## Treatmnt2.5 -0.568 0.486 0.486
## Treatment3 -0.811 0.974 0.974 0.486
## Treatment4 -0.810 0.973 0.973 0.486 0.973
## Temperatr27 -0.179 0.001 0.001 0.000 0.001 0.001
## Temperatr30 -0.174 -0.001 -0.001 0.000 -0.001 0.000 0.453
plot(allEffects(rel_hsat_model_hypnea))
tab_model(rel_hsat_model_hypnea, show.intercept = TRUE, show.se = TRUE, show.stat = TRUE, show.df = TRUE, show.zeroinf = TRUE)
| supersat_rel | ||||||
|---|---|---|---|---|---|---|
| Predictors | Estimates | std. Error | CI | Statistic | p | df |
| (Intercept) | 65.01 | 3.73 | 57.66 – 72.36 | 17.42 | <0.001 | 274.00 |
| Treatment [1] | -7.66 | 4.36 | -16.24 – 0.93 | -1.76 | 0.080 | 274.00 |
| Treatment [2] | -5.66 | 4.36 | -14.25 – 2.92 | -1.30 | 0.195 | 274.00 |
| Treatment [2.5] | -10.80 | 6.22 | -23.06 – 1.45 | -1.74 | 0.084 | 274.00 |
| Treatment [3] | -7.92 | 4.36 | -16.51 – 0.66 | -1.82 | 0.070 | 274.00 |
| Treatment [4] | -10.50 | 4.36 | -19.09 – -1.91 | -2.41 | 0.017 | 274.00 |
| Temperature [27] | 4.54 | 1.41 | 1.77 – 7.32 | 3.22 | 0.001 | 274.00 |
| Temperature [30] | 6.12 | 1.31 | 3.54 – 8.70 | 4.67 | <0.001 | 274.00 |
| Random Effects | ||||||
| σ2 | 24.18 | |||||
| τ00 Plant.ID | 15.73 | |||||
| τ00 Run | 24.72 | |||||
| τ00 RLC.Order | 0.85 | |||||
| ICC | 0.63 | |||||
| N Run | 8 | |||||
| N Plant.ID | 144 | |||||
| N RLC.Order | 6 | |||||
| Observations | 286 | |||||
| Marginal R2 / Conditional R2 | 0.235 / 0.717 | |||||
Construct null model to perform likelihood ratio test REML must be FALSE
hypnea_relhsat_treatment_null <- lmer(formula = supersat_rel ~ Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = hypnea, REML = FALSE)
hypnea_relhsat_model2 <- lmer(formula = supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = hypnea, REML = FALSE)
anova(hypnea_relhsat_treatment_null, hypnea_relhsat_model2)
## Data: hypnea
## Models:
## hypnea_relhsat_treatment_null: supersat_rel ~ Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order)
## hypnea_relhsat_model2: supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order)
## npar AIC BIC logLik deviance Chisq Df
## hypnea_relhsat_treatment_null 7 1887.1 1912.7 -936.53 1873.1
## hypnea_relhsat_model2 12 1870.3 1914.2 -923.17 1846.3 26.723 5
## Pr(>Chisq)
## hypnea_relhsat_treatment_null
## hypnea_relhsat_model2 6.458e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
hypnea_relhsat_temperature_null <- lmer(formula = supersat_rel ~ Treatment + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = hypnea, REML = FALSE)
hypnea_relhsat_model3 <- lmer(formula = supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order), data = hypnea, REML = FALSE)
anova(hypnea_relhsat_temperature_null, hypnea_relhsat_model3)
## Data: hypnea
## Models:
## hypnea_relhsat_temperature_null: supersat_rel ~ Treatment + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order)
## hypnea_relhsat_model3: supersat_rel ~ Treatment + Temperature + (1 | Run) + (1 | Plant.ID) + (1 | RLC.Order)
## npar AIC BIC logLik deviance Chisq Df
## hypnea_relhsat_temperature_null 10 1886.6 1923.1 -933.28 1866.6
## hypnea_relhsat_model3 12 1870.3 1914.2 -923.17 1846.3 20.213 2
## Pr(>Chisq)
## hypnea_relhsat_temperature_null
## hypnea_relhsat_model3 4.081e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Plots
hypnea %>% ggplot(aes(treatment_graph, supersat_rel)) +
geom_boxplot(size=0.5) +
geom_point(alpha = 0.75, size = 3, aes(color = Temperature), show.legend = TRUE) +
labs(x="Treatment (salinty/nitrate)", y= "Relative Hsat (%)", title= "D", subtitle = "Hypnea musciformis") +
scale_x_discrete(labels = c("35ppt/0.5umolN", "35ppt/14umolN", "28ppt/27umolN", "28ppt/53umolN", "18ppt/53umolN", "11ppt/80umolN")) +
ylim(20, 100) + stat_mean() +
geom_hline(yintercept=60, color = "red", size = 0.5, alpha = 0.5) +
scale_color_manual(values = c("#9C0627", "#BB589F", "#F4B4E2")) +
theme_bw() +
theme(legend.position = c(0.88,0.88), plot.title = element_text(face = "bold", vjust = -15, hjust = 0.05),
plot.subtitle = element_text(face = "italic", size = 14, vjust = -18, hjust = 0.05))
For regression with growth Hypnea remove hm6-4 on 11/12 that had no d9 RLC (final weight 0.1017) and hm6-4 on 10/29/21 because it was white and also looked dead
gr_hypnea <- subset(growth_rate, Species == "Hm" & final.weight != 0.1017 & growth_rate_percent > -87.96837)
hypnea$growth_rate <- round((gr_hypnea$final.weight - gr_hypnea$Initial.weight) / gr_hypnea$Initial.weight * 100, digits = 2)
Plot a regression between the photosynthetic independent variables of interest and growth rate
hypnea_growth_rel_hsat_graph <- ggplot(hypnea, aes(x=supersat_rel, y=growth_rate)) +
geom_point(alpha = 0.5, size = 3, aes(color = Treatment), show.legend = TRUE) +
geom_smooth(method = "lm", col = "black") + theme_bw() +
labs(title = "B", subtitle = "Hypnea musciformis", x = "Mean Rel. Hsat (%)", y = "growth rate (%)") +
stat_regline_equation(label.x = 65, label.y = 200) +
stat_cor(label.x = 65, label.y = 215) +
theme(plot.title = element_text(face = "bold", vjust = -15, hjust = 0.05),
plot.subtitle = element_text(face = "italic", vjust = -20, hjust = 0.05))
hypnea_growth_rel_hsat_graph
## `geom_smooth()` using formula = 'y ~ x'
Summarize the means for relative Hsat
hypnea %>% group_by(Treatment) %>% summarise_at(vars(supersat_rel), list(mean = mean))
## # A tibble: 6 × 2
## Treatment mean
## <fct> <dbl>
## 1 0 68.6
## 2 1 61.7
## 3 2 63.6
## 4 2.5 57.8
## 5 3 61.4
## 6 4 58.7
hypnea %>% group_by(Run) %>% summarise_at(vars(supersat_rel), list(mean = mean))
## # A tibble: 8 × 2
## Run mean
## <fct> <dbl>
## 1 1 64.5
## 2 2 67.8
## 3 3 60.2
## 4 3.5 54.2
## 5 4 56.0
## 6 7 57.8
## 7 8 67.0
## 8 9 70.1
hypnea %>% group_by(Run) %>% summarise_at(vars(day_length_avg), list(mean = mean))
## # A tibble: 8 × 2
## Run mean
## <fct> <dbl>
## 1 1 654.
## 2 2 645.
## 3 3 626.
## 4 3.5 621.
## 5 4 611.
## 6 7 688.
## 7 8 635.
## 8 9 636.